Square billiard with a magnetic flux 



O 

o 
o 



Q 

U 
d 



> 
ON 

o 
o 
m 
O 
o 
o 



> 
S3 



R. Narevich 1,2 , R. E. Prange 1 , and Oleg Zaitsev 1 
1 Department of Physics, University of Maryland, College Park, Maryland 20742 
2 Department of Physics and Astronomy, University of Kentucky, Lexington, Kentucky 40506 

Eigenstates and energy levels of a square quantum billiard in a magnetic field, or with an Aharonov- 
Bohm flux line, are found in quasiclassical approximation, that is, for high enough energy. Explicit 
formulas for the energy levels and wavefunctions are found. A number of interesting states are 
shown, together with their wavefunctions. Some states are diamagnetic, others paramagnetic, still 
others both dia- and paramagnetic. Some states are strongly localized. Related systems and possible 
experiments are briefly mentioned. 



PACS number(s): 03.65.Sq, 03.65.Ge 
I. INTRODUCTION 

The behavior of mobile charged particles confined to 
some region and also subjected to a magnetic flux has 
interested physicists since the discovery of the Hall ef- 
fect over a century ago. Quantum effects turned out to 
be subtle and surprising, as attested by Landau diamag- 
netism and the Aharonov-Bohm effect jij . About twenty 
years ago, with the advent of the quantum Hall effect 
and mesoscopic systems, the two dimensional case be- 
came prominent. At about the same time, the develop- 
ment of the subject of "quantum chaos" also focussed 
interest on such systems as being among the simplest of 
the "Gaussian unitary ensemble" universality class. 

There is thus a long history of work on the confined 
quantum motion of charged particles in a magnetic flux. 
It is remarkable that all previous workers overlooked the 
fact that many fundamental cases, some of which have 
been extensively studied numerically, can be solved and 
classified analytically to good approximation. Moreover, 
the results have an interesting and suggestive complexity. 

In this paper we obtain good approximate solutions to 
a couple of simply posed and well studied problems of this 
type. A preliminary version has appeared electronically 
g. Rather than stress the generality of our method, we 
focus on a typical problem: a two dimensional charged 
particle confined to a square billiard in a perpendicular 
magnetic flux. Certain conditions on the flux are required 
to justify the approximations, and we also require the 
energy of the particle to be large. These approximate 
solutions are compared to numerical solutions. The two 
flux configurations considered are a uniform flux and an 
Aharonov-Bohm flux line. It is crucial that the square 
billiard is integrable. Systems other than the square to 
which our methods apply will be mentioned at the end. 

There are many recent research papers in which the 
basic system studied is a square or rectangular billiard 
with a magnetic flux. A number of these are inspired by 
the experiments of Levy et al. 0,0] , which measure the 
magnetic susceptibility of a collection of a considerable 
number of mesoscopic, two dimensional metallic systems, 
each approximately a square. The field is nearly uniform 
over the square in this case. 

In the presence of a magnetic flux, there is the pos- 
sibility that persistent currents exist. In other words, 
it is possible that the equilibrium state has a nontrivial 



current. Indeed, we find that eigenstates of the quan- 
tum system have interesting current densities. Although 
the wavefunction has a fairly simple representation, the 
currents can be quite complex. Some states are predom- 
inantly paramagnetic, others are predominantly diamag- 
netic. Still others may support both strong paramagnetic 
and strong diamagnetic currents which in total nearly 
cancel. However, these states strongly affected by the 
field are rather rare and most states have weak persis- 
tent currents. 

Ideas from the field of quantum chaos have also mo- 
tivated much work. Since the magnetic field intuitively 
has circular symmetry, which 'conflicts' with the symme- 
try of the square, one might expect chaos to ensue ||, 
as is the case with the Sinai billiard. Another theme of 
quantum chaos is that of energy level and wavefunction 
statistics ||. These statistics depend on whether time 
reversal symmetry (or other antiunitary symmetry) is in 
force. A natural way to break time reversal is by a mag- 
netic flux. [The square with a uniform flux still has an 
antiunitary symmetry, however.] 

Diffraction effects, in which a classical length shorter 
than the wavelength becomes important, are much stud- 
ied in this context. This is obviously the case for the zero 
radius Aharonov-Bohm flux line . The sharp corners of 
the square also cause diffractive effects in the presence of 
a uniform magnetic field. These effects are much smaller 
than for the flux line, of course. We give estimates for 
the parameter range in which such diffraction becomes 
important, although we defer study of these effects. 

Our main motivation however, is that we add to the 
store of solvable problems, and perhaps suggest some ex- 
periments. In the textbooks, there are relatively few such 
integrable problems, basically, only those which reduce to 
one dimension, or separate into several one-dimensional 
problems. The square without a magnetic flux is such a 
case in which the x and y motion separates. 

The traditional way to widen the class of approxi- 
mately solvable problems is perturbation theory and in- 
deed, our approach is a form of quantum perturbation 
theory. We are able to study systems which are inte- 
grable except for a 'classically weak' perturbation. Of 
course, treating weak perturbations classically is chal- 
lenging, because the long time behavior may be chaotic. 
However, quantum perturbation theory is better behaved 
and depends on the short time rather than the long time 
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classical behavior. 

On the other hand, perturbations small in this sense 
can give rise to very large quantum effects, especially on 
the wave functions. Moreover, even perturbations which 
change some length scale by an amount 5L << A, where 
A is the wavelength, can have big effects. In terms of the 
standard nomenclature, a degenerate perturbation the- 
ory is required, and a number of unperturbed states are 
strongly mixed together to give the final result. Usually, 
this is done by diagonalization of a small matrix, but in 
our case the 'matrix' can be quite large. However, rather 
than just diagonalizing some matrix by a computer cal- 
culation, we obtain the result by an intuitively appealing 
Schrodinger differential equation. 

Having solved, for the first time, this interesting class 
of problems, we discovered some other methods of ob- 
taining the solution at the same level of approximation. 
We shall present these methods elsewhere S. 



II. SQUARE IN A UNIFORM FIELD 

We begin with the case of uniform field. For given ve- 
locity v, the cyclotron radius is R c = v/uj c = cp/eB 
where oj c — eB/mc, and p — mv. The momentum 
p = hk = h/X is quantally related to wavenumber k 
and wavelength A. We define the classical small param- 
eter e = L/R c = eBL/hck = 2ircf>/(f>okL. Here L is the 
length of the side of the square, <j) is the magnetic flux 
BL 2 and <fio is the flux quantum hc/e. Small e allows us 
to approximate orbits within the square as straight lines, 
to first approximation. This is sometimes known as the 
Aharonov-Bohm regime ||, since the leading quantum 
effects come from the phase interference effects associ- 
ated with the vector potential, and do not depend on the 
change of classical orbital motion caused by the Lorentz 
force. Many potential experiments are in this parameter 
range. 

We choose units such that the dimensionless field is 
B = 2"K(f)/4>Q, i-e., 2tt times the number of flux quanta in 
the square. We take L, fi and 2m to be unity so that 



e = B/k « 1. 



(1) 



The dimensionless wavenumber k is the number of wave- 
lengths in a side of the square, up to a factor of 2ir. It 
satisfies 



k » 1, 



(2) 



which is the basis for the quasiclassical approximation. 

We shall see that the condition for standard quantum 
perturbation theory to work is 



ky/e « 1, 



(3) 



or, in other words, ykB << 1. This is not completely 
obvious, and in other contexts |J, it has been guessed 
incorrectly that the quantum 'perturbation border' is 
ke << 1, i.e. B « 1, which has the simple meaning 
that the number of flux quanta in the square is small. 



We, however, find that nothing much changes at the bor- 
der B ~ 1. 

There may also be a 'high energy' condition in the form 
of a requirement that ke b << 1. The exponent b depends 
on the smoothness of the perturbation. We find that for 
the uniform field, 6 = 2, while for the ideal flux line b = 1. 
Note that for fixed B the energy can be arbitrarily high, 
but if instead e is kept fixed, there is a limitation on the 
energy. 

This high energy condition is basically the requirement 
that diffraction effects not be too important. From a 
semiclassical perspective, diffraction effects occur where 
the classical system has a length scale as short or shorter 
than A. An ideal flux line obviously gives rises to diffrac- 
tion effects. Billiards also have such a length scale of 
course, namely the distance it takes the confining poten- 
tial to change from zero inside the billiard, to infinity 
outside the billiard. This can be taken into account by a 
'Maslov phase' of w at the boundary, however. There are 
also the sharp corners of the square. The square corners 
whose angular opening is n/N, where N = 2 is an inte- 
ger, are a special case at which no diffraction occurs [[to). 
With such an angle, the billiard can be extended by re- 
flection, and the corner in effect disappears. However, in 
the presence of a magnetic field, this reflection technique 
does not work, and with sufficiently large field, orbits 
which hit directly into the corner eventually become im- 
portant to the semiclassics. 



III. BOGOMOLNY'S QUASICLASSICAL 
SURFACE OF SECTION METHOD 

Our approach Jll| utilizes the quasiclassical surface of 
section [SS] method of Bogomolny |12|. Poincare's sur- 
face of section is a surface in classical phase space through 
which all interesting orbits repeatedly pass. For two di- 
mensional systems, the surface of section is a two dimen- 
sional phase space. For a billiard, the Birkhoff surface 
of section is often chosen. Namely, the space part of the 
surface of section represents a point on the boundary at 
which the orbit bounces and is usually measured by the 
distance along the boundary of the billiard. The variable 
conjugate to this is the component of momentum paral- 
lel to the boundary at the moment of contact. However, 
many possible surfaces of section can be considered, and 
some are more convenient than others. 

Bogolmony's method is a generalization of the "bound- 
ary integral method" |T^,[l3|], applicable for billiards, and 
based on Birkhoff 's surface of section, to much more gen- 
eral systems and surfaces of section. The boundary in- 



tegral method introduces an operator K(x, x' , E) and an 
integral equation tp (x) — J dx'K(x, x'\ E)tp(x'). This ex- 
act equation has nontrivial solutions only when E is on 
the spectrum. The SS wavefunction tp(x) is the normal 
derivative of the full wavefunction, ip(x) = d^>(r)/dn, 
when r is at the boundary point x. We should mention 
that only recently has the boundary integral method been 
extended to uniform magnetic fields in the case that e is 
of order unity fla] . 
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Bogolmony's operator T(x, x'; E) is basically the qua- 
siclassical approximation to K. It thus takes the particle 
crossing the SS at position x' to its next crossing at posi- 
tion x, all at energy E = k 2 . The quasiclassical approxi- 
mation to the spectrum is determined by the existence of 
solutions of Tip — ip. If only the spectrum is of interest, 
as it has been for many authors, the condition may be 
expressed as det [1 - T(E)] = 0. 

The operator T is given quite generally by 



T(x, x; E) 



1 

2tt7 



d 2 S(x,x']E) 



dxdx 1 



exp [iS(x, x ; E)j 



(4) 



where S = f , p • dr is the action integral along the clas- 
sical path from x' to x. Note that, rather than giving 
position and momentum on the SS, positions at two se- 
quential crossings of the SS are given. It is assumed, for 
notational convenience, that there is a unique orbit from 
x' to x. We also suppress the Maslov phase. Note that T 
is semiclassically unitary. 

The classical action S{x,x') generates the surface of 
section map. Namely, the momenta conjugate to x, x' 
are given by 



dS{x,x') 
dx 



dS{x,x') 
dx' ' 



(5) 



Eq. (|^) implicitly gives the surface of section map 
{p, x}=M{p',x'}. 



of width 2 obtained by reflecting the original square first 
about x = i and then about y — ^ , and finally repeating 
the resulting 2x2 square periodically to x — ±oo. The 
flux changes sign in neighboring squares. This geometry 
is shown in Fig. 1. There are a continuum of channel so- 
lutions. The solutions to the original square are a subset 
of these, which exist only at certain quantized energies. 
This quantization can be carried out in several ways, one 
of which is shown below. 

The SS is taken as the axis y — — \ which is identified 
with y = |. Because the field is classically weak, the path 
used to calculate the action is approximated by a straight 
line. We immediately find 



S(x, x') = k\jA+{x-x'f + $(x, x'), 



(6) 



i. e. the flux free result plus $ = (e/c) J A-dr, where the 
integral is done along the straight line path. 

Our scheme jll| finds solutions of Tip = ip by a per- 
turbation theory. Taking advantage of the fact that T is 
unitary, we first solve 



dx'T(x,x';E)ip(x') = e iu ^ip{x), 



(7) 



treating k = \f~E~ as a parameter, and then find the en- 
ergies by solving uj(k) = 2nn. Given ip(x), a quadrature 
which can be carried out quasiclassically yields the full 
wave function ^(x, y). Details are given in Appendix A. 



v = 1/2 




FIG. 1. A channel of 2x2 squares replacing the original 
lxl square. Adjacent squares are reflected and the magnetic 
flux changes sign. An orbit in the original square is replaced 
by a nearly straight line orbit in the channel. The orbit cur- 
vature of radius R c due to the Lorentz force is exaggerated 
for clarity in the figure. The orbit shown goes from x' to x + 2 
in the channel representation and is close to a (1, 1) periodic 
orbit. It is paramagnetic, since the reflections from the square 
sides cause it to circulate in the opposite direction from a free 
orbit in the field. 

We attempt to simplify T by astute choice of the sur- 
face of section. It seems simpler to use just one side of 
the square, rather than all four sides. It is even easier to 
use a method of images. Namely, we consider, instead of 



IV. RESONANCES 

Classical and quantum perturbation expansions in 
powers of e fail near resonances, necessitating modifi- 
cations which introduce y/e. Classical resonances cor- 
respond to periodic orbits of the unperturbed system. 
Periodic orbits on the square correspond to straight line 



orbits in the channel from (a/, — i) to (x 



■2p/q, §). 



nary perturbation theory works. See Ref. fill] for a fuller 



p a : 

m 



a unit square, x, y £ 



] , an infinite channel 



Here q is a positive integer and p is a positive or nega- 
tive integer relatively prime to q. Negative and positive 
p are not equivalent if there is a magnetic flux. Between 
resonances, or near resonances with large p and q, ordi- 
nary pertui 
discussion. 

We now specialize to the (±1, 1) resonances. These are 
the simplest resonances depending strongly on the field 
and, as we shall see, in some sense dominate the magnetic 
response. 

We look for a solution of Eq. (|7]) of the form tp(x) = 
e lKX u m (x) where k — fccos45° = k/y/2 and u m varies 
much more slowly than the exponential. The reason for 
this choice is that the phase factor e lKX makes the rapidly 
varying phases in the integral J dx'Tip to be stationary 
at x' — x — 2. This corresponds to rectangular shaped 
periodic orbits of the original square whose sides make 
angles of 45° with the x axis. Such orbits are shown in 
Fig. 1. 

Because e is small, the phase $ does not greatly influ- 
ence the position of the stationary phase and it suffices 
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PI to evaluate <§>(x,x') at $(x,x — 2) = $(x + 2,x). 
[The accuracy of this approximation depends on the 
smoothness of $, and failure of the approximation is re- 
lated to the onset of diffraction effects mentioned earlier.] 
<I>(a; + 2, x) is obtained by integrating the vector potential 
about the closed rectangular loop, and somewhat remark- 
ably is independent of gauge. The result |yj is that u m 
satisfies the Schrodinger equation 

- u"> + V(x)u m = E m u m (8) 

where V(x) — — fc$(x + 2, x)/C and C = v° is the length 
of the periodic orbits. Thus we convert the phase $ to a 
'potential' V. The transverse energy E m enters into the 
total energy of the eigenstate, according to Eq. (|T^) be- 
low, and m is one of two quantum numbers classifying 
the states. 



V. TRANSVERSE POTENTIAL AND PERIODIC 
ORBITS 

For a given resonance and surface of section, there 
is an effective potential V which determines the func- 
tions u m and the energy E rn . The resonance classically 
corresponds to a continuous set of nonisolated periodic 
orbits of the integrable problem. Before perturbation, 
each of these orbits has the same action. The poten- 
tial V, to leading order, is proportional to the change 
of the action under perturbation, calculated along the 
unperturbed path. Each such path is labelled by the pa- 
rameter x, where it crosses the surface of section. Higher 
order corrections may also be found |fl[] . 

Knowledge of the potential gives much qualitative in- 
sight into the problem. Its minima, [if smooth], are at 
stable periodic orbits, as a rule, and its maxima are at 
the unstable orbits. In that sense, it represents a classical 
island chain. Of course, it was known how to quantize 
states near the stable periodic orbits, if a harmonic ex- 
pansion is allowed. However, the states that can be found 
with the aid of V are much more general and in partic- 
ular the states with energies E m near or even above the 
maxima of the potential can also be found. 

In general, isolated unstable periodic orbits do not 
support wavefunctions, but rather 'scar' them | Ejj |. In 
other words, there appears some excess weight on the 
wavefunction near the unstable orbit. In the sense 
of Feynman's path integral formulation, there are not 
enough classical paths 'near' the unstable orbit, to build 
a complete wavefunction. Here 'near' means that the 
paths are close to the periodic orbit in the sense of being 
well approximated by a quadratic expansion about the 
periodic orbit. 

The same is true in the present case, and wavefunc- 
tions cannot be built just from orbits near an unstable 
periodic orbit. However, because of the small parameter 
e, we can approximate well an entire shell of orbits in 
the Feynman integral, and express the result in terms of 
the potential V(x). This shell can support many states 
which we find. There are states whose energies are near 



the maxima of V and thus have extra weight near the 
unstable periodic orbits. 

The interpretation of u m is that it gives the structure 
of the wavefunction 'transverse' to the resonant periodic 
orbits. Along the periodic orbits, the wavefunction varies 
rapidly, but transversely, it varies relatively slowly. The 
'longitudinal' and transverse motions are weakly coupled, 
because V and thus u m and E rn depend on fc, but this is 
easy to take into account. 

The concept of 'transverse' is a little murky in the 
quantum case, although there are cases, including the 
one under study where it can be made more precise. We 
shall not dwell on this further in this paper, however. 

We also remark that the Schrodinger equation (||) re- 
quires boundary conditions, in order to pick out the phys- 
ically interesting solutions. These boundary conditions 
come from the properties imposed on the solution by the 
physics of the problem, and are usually simplified by sym- 
metries of the problem. 



VI. UNIFORM FIELD SOLUTION 
A. Effective potential 

For the uniform field, the potential is 

V(x) = -Bk (| - 2a; 2 ) /£; x€[-±,±], 

V(x) = +Bk [i - 2{x + if] /£; x E [-§,-£], 

V(x) = V(x + 2). (9) 

The factor | — 2a; 2 is simply the area enclosed by a pe- 
riodic (1,1) resonant orbit in the shape of a rectangle 
which bounces from the bottom of the square at x. In 
Appendix B we give the corresponding potential for other 
resonances. 

This periodic potential consists of alternating positive 
and negative harmonic potential wells of depth Bk/2£. 
At the boundaries x = the second derivative of the 
potential is discontinuous, a fact which leads to the men- 
tioned diffraction effects at sufficiently large B 2 /k. 

For the (—1,1) resonance, whose orbits are time re- 
versed (1, 1) orbits, V(x) changes sign. This would not 
be true if V had its origin in a time reversal invariant 
perturbation of the square, for example, a small change 
of shape. We can include the (—1,1) resonance in the 
present scheme by attributing the region 1/2 < x < 3/2 
to that resonance. This extension of the x coordinate is 
thus similar to use of a 'angle' variable, with positive x- 
velocity v x for x £ [— i, i], and negative v x for x € [h, §]. 

If v 7 Bk = k^JT is small, the potential V(x) can be 
treated perturbatively. On the other hand, for suffi- 
ciently large Bk/C, Eq. (||) will have low lying tight 
binding harmonic oscillator type solutions centered at 
x = 0, (if B > 0),with energies approximately given 
by 

E m = -\Bk/C + {m+ ±)y/8Bk/£. (10) 
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This formula holds for m < < W Bk/C. The lowest wave- 
function is approximately uq(x) = e ~V Bk / 2Cx which is 
arbitrarily narrow at large energy. These states are para- 
magnetic, as follows from the fact that dE m /dB < 0. 
This will be seen more clearly below. 

Eq. (@) is valid for larger m. Although very simple an- 
alytic answers are not available, the problem is the well 
known one of a particle in a one dimensional periodic 
potential. We shall see below that we need only con- 
sider the boundary conditions u(x + 2) = ±u(x). This 
simplification is a consequence of the symmetry of the 
square, and something slightly more complicated would 
be needed for the rectangle. 

The solutions to Eq. m) may be put into four classes, 
A, B,C, D. Class A states are those with 'low' energies 
near the bottom of the well, E m ~ — ^Bk/C For these 
cases u(x) has support only near x = 0, ±2, ±4, ... These 
localized states are strongly paramagnetic, that is, the 
current circulates in the opposite direction from that of 
the particle in free space. In this case dE m /dB < 0. [We 
shall see that the transverse energy E m carries nearly 
all the field dependence of the total energy of the corre- 
sponding two dimensional eigenstates.] 

Class D states have energies much greater than the 
maximum potential energy, i.e. E m » \Bkj L. In this 
case, the magnetic field is a small perturbation, since 
the 'potential energy' V(x) in Eq. (B) is small compared 
with the 'kinetic energy' given by u". These states are 
weakly diamagnetic. We shall not consider further this 
case. Of course, the approximation of expanding about 
the (1, 1) resonance eventually breaks down, and higher 
order resonances are eventually involved 1 1 1 . 

Class C has total transverse energy near the top of 
the potential V(±l), that is, E m ~ ^Bk/C and the 
states are strongly affected by the magnetic field. Very 
crudely, they are somewhat localized or 'scarred' near 
x = ±1, since they spend more time in that region. This 
means that they are strongly influenced by the (—1, 1) 
resonance. They are diamagnetic and dE m /dB > 0. 

States of class B form a transition region between 
the low-energy paramagnetic states, and the higher en- 
ergy diamagnetic ones, i.e. near E rn ~ 0. These are 
states strongly affected by the field, but are such that 
dE m /dB ~ 0.' 



B. Quantization 

There are two states with identical energy in the re- 
peated square scheme. These are ipi — e tKX u m (x) 7 and 
ipu = e~ lKx u m {x — 1). [Changing the sign of n is equiv- 
alent to changing the sign of the field, which in turn can 
be accomplished by replacing V(x) by V(x + 1).] 

Rather than finding the eigenvalues by imposing con- 
ditions directly on the ^'s, as in the Appendix A, we pro- 
duce the four two-dimensional solutions , 3> corresponding 
to 4>i,ii- [Each ip{x) gives two \l/(x,2/)'s because y and 
1 — y in the strip represent the same point in the origi- 
nal square.] We show elsewhere || that these states can 



also be found directly by a Born-Oppenheimer approxi- 
mation. One of these states may be written 



* (x,y) = e iK ^u m (x-y-±) 



(11) 



The remaining three states, 1,2,3, can be obtained by 
rotations, e.g. *5!\(x,y) = TZ 1 i'o(x,y) = *f?o(y, — x), etc. 
Here TZ : (x,y) — > (y, — x) is the rotation by 90°. The 
gauge can be chosen so that the Hamiltonian is invariant 
under TZ. Therefore, the symmetry of an eigenstate can 
be labelled by r — 0, 1, 2, 3, where the eigenvalue of TZ is 
i~ r . 

Thus, an eigenfunction with symmetry r is given by 



{x-y-h). (12) 



[The sequence of rapidly varying phase factors is 
| e i«(a+ ! /) ) e «(-x+y) ) e - lK (x+y)^ e M^-v)j _ These in turn 

are rapidly varying in the 45° directions of the sides of the 
periodic orbits.] In general, a solution of Eq. (||) satisfies 
the boundary condition u m (x + 2) = e l ^u m {x). We need 
to find the allowed values for /3 and n which will give the 
quantized energies. These conditions are obtained by re- 
quiring ^r r \(x, — ^) to vanish, corresponding to Dirichlet 
conditions in the original problem. If the wavefunction 
vanishes on the bottom, it will by symmetry vanish on 
the boundary of the square. 

Clearly, the sum of the two terms [s = 0,3] in 
Eq. (jig) which are proportional to e +ltiX must van- 



ish. This implies u m {—x) = —i e 



t (x). The re- 



flection symmetry V{x) = V(—x) allows us to take 
Um(x) = (—l) m u m (—x). In turn, this allows quantiza- 
tion of k in the form k = nir/2, where n is an integer satis- 
fying certain conditions depending on r and m. Similarly 
the two terms proportional to e~ tKX in Eq. ( |l2| ) give the 
condition e 1 ^ — (— l) r . The relationship of n to r and m 
is 

nmod4 = [2(1 - mmod2) + r] mod 4. (13) 
It is straightforward to find for the eigenwavenumber 



k n 



2irn/£+E m /k. 



(14) 



Eq. (jl4|) should be solved iteratively. For example, the 
first approximation replaces the k dependence of the term 
E m /k by 2itn/C Equivalently, the energy 



E n , m = An 2 n 2 /C 2 



2E„ 



(15) 



Note that E m depends, relatively weakly, on n, since 
the k in formula Eq.(jlO|) should be replaced by 2-kti/C. 
The dependence of the total energy on B comes through 
the term E rn . Eqs. ( |l4| ) and (jl^) hold for all symmetries 
and successive values of n at fixed m cycle through the 
representations of 1Z. Note that, since E m /k << k the 
wavelength is given approximately by C/n, i.e. the length 
of the classical orbits is an integer number of wavelengths. 

Thus we have an expression for the energies of a class 
of states, namely the (±1,1) resonant states. They are 



5 



labelled by integer n which effectively give the number 
of wavelengths measured along the (1,1) periodic orbits, 
and by a second integer m which gives the number of 
'nodes' 'perpendicular' to this orbit. The very low m 
states could very well have been found by earlier meth- 
ods, since they can be obtained by expansions about the 
stable periodic orbits. However, these remarkable states 
do not seem to have been noticed heretofore. 



C. Orbital magnetism 

The (1, 1) states just obtained dominate the magnetic 
orbital susceptibility in a parameter range appropriate 
to experiments ||. The susceptibility for the square is 
on a scale rather larger than the Landau diamagnetism. 
It is of course not necessary to find the states, or for 
that matter, their energies, to calculate the susceptibil- 
ity. That is because the susceptibility depends only on 
the density of states smoothed over an energy width pro- 
portional to the temperature. The Gutzwiller or better, 
the perturbed Berry- Tabor trace formula Q] is designed 
to give exactly that quantity in quasiclassical approxima- 
tion. Nevertheless, it's interesting and previously unre- 
marked, that a small subset of states accounts for most 
of the magnetism. 

We start by finding the orbital susceptibility \ of a 
system of nonintcracting electrons in a grand canonical 
ensemble. This is given by \ = / dB where the mag- 
netization M. = —dfl(T,/i,B)/dB. Here the grand po- 
tential is 



n(T,/i,B) = -fc B T^ln \\ + e -^-^/ kBT 



(16) 



The temperature is T, k B is Boltzmann's constant, and 
H = kp is the chemical potential. The dependence of £1 
on B comes only because the eigenenergies E a depend on 
B. The sum is over all eigenstates labelled by a. 

We divide the states a into those relatively few whose 
energies depend appreciably on the field and the rest. 
These field dependent states are exactly the (1,1) states 
found above, plus possibly states classically associated 
with a few other low resonances, e.g. (1,3). The reason 
for this is that the (1,1) states enclose the maximum 
directed area. They also have the shortest length C which 
we will see plays a role. The even shorter (0,1) periodic 
orbits do not enclose any flux in the approximation of 
neglecting the curvature of the orbits although at higher 
fields they eventually become important. 

Thus, we replace J2 a = St + J2 n ,m anc ^ we can nc ~ 
gleet the sum b over field independent states. The second 
sum, over (1, 1) resonance states, has many fewer terms 
than the first in a given range of energy. Since [i is re- 
lated to the number of particles, it is nearly independent 
of B. It is possible to find \i = k F + 5/j,(B) and make a 
consistent expansion, and that is indeed necessary if an 
average over a large number of squares with canonical 
statistics is done Q. However, we just want to illustrate 
how the (1,1) states dominate the susceptibility, and we 



will not consider this further average. Then, we may 
approximate 

f)E 

M(T,»,B) = -^-£pf D [E n , m (B)]. (17) 



and fo is the Fermi-Dirac distribution function. 

Using the Poisson sum formula, replace the sum on n 
in Eq. ( |l7j ) by an integral over k, and do the integral to 
obtain 



M 



k B T 
kF 



^ a m C exp 



2k f 



x sin 



Cs ( kF 



E,, 



(18) 



Here, ui r = ir(2r + l)k B T is the Matsubara frequency and 
a m = dE m /dB. [We have dropped the 'leading' term in 
the Poisson formula which totally neglects the discrete 
quantum nature of the states and which therefore cannot 
produce a magnetization.] As an example ||f|], take k B T 
ten times the level spacing of all levels, i.e. k B T = 2Qir 
in our units. Then, uj C/2 « 300. If k F « 300 - 600, 
so that the square contains about 2-6xl0 4 electrons, the 
exponential suppression will not be too serious for r — 
0, s = 1. However, larger r or s do not contribute much. 
[In the trace formula approach, s gives the number of 
repetitions of the primitive periodic orbit and the sum 
over r is explicitly carried out.] 

Eq. (|l^) shows that states with larger £, i.e. smaller 
spacing, are suppressed, exactly as seen from the trace 
formula in terms of periodic orbits. It also shows that 
relatively large field dependence of the levels, a m , is im- 
portant. For the square, the (1,1) states have the small- 
est C and also the largest a m . The (2, 1) resonance does 
not couple to a small constant field. 

It is also seen that if Eq. (|l^) is averaged over many 
squares of somewhat different sizes, because of the oscil- 
lations of the sine, the result is much reduced and it is 
necessary to go to higher order in 8fi. 




-0.5 - 



FIG. 2. 8E m /dB vs E m /B. Stars are for 

m = 0, 2, . . . m m ax = 14, . . . and n = 62, B = 25. The contin- 
uous curve is Eq. (^) Q f the text. 
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In Fig. 2 we show dE m /dB as a function of E m /B 
for the (1,1) resonant states, where B = Bk/2£. Ac- 
cording to Eq. (|l8|), if the sign of dE m /dB is nega- 
tive the contribution of the corresponding state is para- 
magnetic. For an exact, two dimensional wavefunction 
^n,m with energy E n _ rn , it is known that dE n ^ m /dB 
= —(£/k) J d 2 r[r x j(x,y)] z , i.e. the expectation value 
of the z-componcnt of the magnetization density. Here 
the current is 



j{x,y) = 2Re^*{x,y) 0V-A(x,y)j 9{x,y). (19) 

In our approximation, according to Eq. ( |l5| ) above, 
OEn^m/dB is equal to 2dE m /dB. It follows from Eqs. 

(§) and (g) above, that 8E m /dB = (v(x) 
V(x) 



xe[- 



where V = V/B. Since V < for 

and V > for x £ [|, |] , etc. we see that 
the sign of the magnetic response of a given wavefunction 
depends on which region of V dominates the expectation 
value. Of course, the classical periodic orbits in these 
two regions have the expected sense. 

Finally we may express dE m /dB quasiclassically as 



dE m 
dB 



J dxV(x) E - V(x) 



Jdx 



E - V{x) 



(20) 



D. Visual representations of eigenstates 

Eigenstates of Hamiltonians not satisfying time rever- 
sal invariance are rarely shown graphically in the litera- 
ture, except for trivial cases. In the invariant case, nu- 
merically produced picture galleries of eigenstates for sys- 
tems such as the Bunimovich stadium p"5| ], have led to 
a great deal of interest and insight, both theoretical H] 
and experimental [ jig ]. 

'Magnetic' states are not real but are inherently com- 
plex. A complete graphical picture of a complex state 
would seem to require twice the number of pictures as 
that necessary for a real state satisfying time reversal in- 
variance. In addition, the states themselves are gauge 
dependent, and only gauge invariant quantities are phys- 
ically meaningful. 

A difference of the magnetic case as compared with the 
time reversal invariant case is that generally \P or j do not 
vanish along nodal lines, but rather only at isolated nodal 
points. There may of course be symmetries or boundary 
conditions requiring these quantities to vanish along a 
line, but in general the representation of wavefunctions 
by their nodal patterns |L5) is not available in the absence 
of time reversal invariance. 

Two gauge invariant quantities we choose to dis- 
play are the absolute value squared of the wavefunction 
y)| and the current density. The current density 
is a two dimensional vector field that is divergence free, 
V • j = 0. The one-dimensional surface of section states 
are also of interest. 



Here E = E m /B can be treated as a continuous variable, 
so that dE m /dB as a function of E falls on a continuous 
curve which in this approximation is independent of B. 
The integrals are between the turning points. For E near 
the minimum V m i n of V , dE m /dB ps Vmin- For E at the 
maximum of V, the integrals diverge at x = 1, thus mak- 
ing dE m /dB positive, and also giving the cusp in Fig. 2. 
Of course, our simple quasiclassical approximation needs 
corrections in this case. 

The quantization of the transverse motion relates E = 
E m /B and m according to the standard formula 



. . \ E - V(x) =7r(m + i)/ 1 



(21) 



The spacing of the quantum states therefore depends on 
B. 



The level mtop such that E m 



V(l) satisfies, ap- 
proximately m top ps 0.6V B. The states with k = 142 
and B — 25 were chosen in part because with these pa- 
rameters there is a state very near the diamagnetic max- 
imum m t0 p = 14, and another state with m = 10 for 
which 8E m /dB is very small. On the other hand, with 
these parameters, while our approximate states are very 
good for small m, the transverse momenta are becom- 
ing large enough (as compared with the longitudinal mo- 
menta) that our approximation deteriorates significantly 
at m > 10, especially in regions of the square where the 
wavefunction is small. 




FIG. 3. A three dimensional view of the absolute value 
squared of localized paramagnetic state for n = 62, B = 31.4 
and m — 0. The cut near the side of the square is proportional 
to uo which is close to a Gaussian. Near the side there are 
interference oscillations from two main terms in Eq. jl2]). 

We first show a picture of a class A state, which is 
quite simple to represent. Fig. 3 shows \^(x, y)\ 2 , for 
n = 62, k fa 2vr62/£ ps 140, B = 31.4, J~Bk ps 66, and 
m = 0, which implies the symmetry r = 0. For such a 
well localized u m , each term in Eq. ( p"2| ) dominates one 
side of the rectangular periodic orbit. For example, near 
x = —y = j only the first term, s — 0, in the sum jl^ ) 
is appreciable. In this region |^(a;, y)\ ps uq(x — y — |) 
which is well approximated by a Gaussian. 
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There is interference near the square edges (e.g. near 
(0, — i)), and two terms contribute appreciably. Near 
this point then 

e^y Uo (x - y - |) + e-^uoi-x - y - |) 
ss |2ito(a;) cos ^niry\ (22) 

so that 1^(0, — 5 + ~)| at its first maximum near (0, — i) 
is about twice as large as |^(|,— [Note that 

n mod 4 = 2.] The current density is thus largest close to 
the middle of the square edges and it is of course nearly 
parallel to the edge there. In this case, the shape of W 
close to an edge is given by u m {x), as can be seen in Fig. 
3. 

Fig. 4 shows the streamlines of the current in this 
state. The direction of the streamline gives the direction 
of the current flow, while the density of streamlines is pro- 
portional to the magnitude of the current density. That 
is, between any two neighboring streamlines, the same 
total current flows. The state in Fig. 4 is paramagnetic, 
that is, the current circulates in the opposite sense from 
that of a free particle in the field. The choice of which 
streamline to display is easily obtained in this case, since 
each line crosses a symmetry line like x — 0, y E [— i, 0], 
once and only once. The current has but one interior 
zero, at the center of the square. 




FIG. 4. Streamlines of the current of the state shown in 
Fig. 3. The interference peaks near x = are apparent. 

Fig. 5. shows the state n = 60, B = 31.4, m = 1, 
which has one 'transverse node'. This node does not give 
rise to a nodal line in the total wavefunction, of course, 
although the wavefunction is small in regions correspond- 
ing to the node. The states have rotational symmetry, 
and we show a different representation in each corner of 
the unit square. 

One remarkable feature of these states, m = 0, 1 is 
that they are localized very near the central paramag- 
netic diamond orbit. Theoretically, for large Bk this lo- 
calization can be as tight as one pleases. Although this 
is a result of our theory which uses a 'potential' function 
V(x), it is clear that V(x) is really something that arises 
from Aharonov-Bohm phases, and there is no classical 
localization based on energy considerations. Although a 



detailed analysis of this can be a little tricky, we find 
the same kind of effect in the AB flux line case, for large 
Bk and small B. This localization is thus like Anderson 
localization jl7]] in the sense that, absent phase interfer- 
ence, localization would not exist. Of course, the random 
disorder aspects of Anderson localization are absent. 




FIG. 5. The state with n = 60, B = 31.4 and m = 1. 
Counterclockwise from lower left, a density plot of 2 , cur- 
rent streamlines, a density plot of |j| , and a dot-stick repre- 
sentation of the current. [The dot is at the calculated point, 
the size and direction of the stick represent the current den- 
sity magnitude and direction.] The function wi(a;) has a node 
at x = 0. This appears as a narrow valley in the center of the 
ridge of current, although there is no node or nodal line there. 
The numerical result is shown, which is very well represented 
by theory. 



V - - ^ ^ — 
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FIG. 6. The state, m = 14, n = 62, B = 25 of maximum 
diamagnetism, as found numerically. The same representa- 
tions as Fig. 5. 

The next figures, Figs. 6-11, have n = 62, B = 25. 
Figs. 6,7 show a state of class C corresponding to m = 14 
which is energetically at the top of the periodic poten- 
tial corresponding to the (—1, 1) resonance. This state is 
diamagnetic with the current circulating in the opposite 
sense from the states with m = 0, 1. The theoretical and 
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numerical wavefunctions are shown. 

The theoretical predictions in this case are consider- 
ably less good than for the states with small m. First, 
the transverse wavenumber is no longer quite so small 
compared with the wavenumber along the path. Second, 
since according to Eq. (fl2|), as many as four approx- 
imately found component states are added, there may 
be relatively large errors, especially in those parts of the 
square where destructive interference is important and 
the final wavefunction is small. Nevertheless, Eq. (|l^) is 
a quite good representation of the more accurate numer- 
ical results. 

In Fig. 8 we show the state m = 10, for which 
dEio/dB is very small. The streamlines of this state 
are rather striking. Note that there are large current 
loops which have opposite magnetic polarity. Again, the 
approximate wavefunction captures many features of the 
exact one, although it docs not reproduce the finer de- 
tails. In this case we show also, in Fig. 9, the transverse 
state uio(x) and the normal derivative of ^62,10 on the 
surface of section. Although for small m, the normal 
derivative on the surface of section bears an understand- 
able relation to u m , it is quite complicated for transverse 
energies this large. 

Finally, in Fig. 10 we show streamlines for the se- 
quence of states m = 6, 8, 10, 12. Although there are sys- 
tematic changes of pattern, we have not tried to ratio- 
nalize these changes. We conclude that even though the 
wavefunction of Eq. ( |l2| ) is fairly simple, it is difficult to 
foresee interference patterns when four terms are impor- 
tant. 




FIG. 7. The same state as Fig. 6. as predicted by our 
approximation. For m this large, the theory does relatively 
poorly, but nevertheless is qualitatively correct. 

These relatively complicated states are harder to rep- 
resent adequately. In both theory and numerics the 
V • j =0 character of the current is not exact. Numeri- 
cally following a streamline, particularly in the neighbor- 
hood of a zero of the current is difficult. We therefore 
imposed the divergence free character by representing 
j(x, y)= V x zx(x, y). We calculated x as a symmetrized 
integral of j x , where j x was obtained either numerically 



or theoretically. The streamlines are then contour lines 
of x- 

The diamagnetic states are not so spatially localized 
as the low m states, although they continue to have a sort 
of localization in 'momentum' space, as we show below. 
More generally, as m increases, the states become more 
delocalized, and eventually become independent of B. 
This means that for larger m all four terms in Eq. ( p^ ) 
make comparable contributions at an arbitrary typical 
point x, y, whereas in the localized case, only one or two 
terms contribute. This gives interference oscillations in 
\ip m \ near |ar| ss \ as shown in Fig. 9. 




FIG. 8. The state, m = 10, n = 62, B = 25 of nearly can- 
celling para- and diamagnetism. Numerical results are shown 
as before, except for the upper left corner where the theo- 
retical streamlines are shown. Theory and numericals differ 
primarily in the low current regions. Notice there is a dia- 
magnetic current loop encircling the diagonals and cancelling 
paramagnetic loops in the triangular wedges between the di- 
agonals. 




FIG. 9. The upper figure shows, for m = 10, 

n — 62, B — 25 the normal derivative \dW 62,10 /dy\ _i , 
obtained numerically. The lower figure is uio(x). Because 
Mio(a;) extends significantly outside the domain [— |, 5] , and 
is 'folded' back into that domain with a rapidly varying addi- 
tional phase factor to construct d^/dn, these two constructs 
are not so intuitively related as for lower m states which are 
more localized. 
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E. Momentum localization 

We have assumed that n >> m, and that u m is slowly 
varying compared with e lvnx / 2 _ We may expand u m = 
X) u m ,ie l7rlx , where / is an integer for r even and half 
odd integer for r odd. Also u m ,i = (— l) m u m ,-i- The 
unperturbed states [sin7rp(x + h)simrq(y + |)] can be 
labelled by integers p, q with unperturbed energies (p 2 + 
q 2 ), dropping a factor it 2 . 

Eq. (|l^) is a superposition of unperturbed states with 
quantum numbers p = + 1, q = \n — I. In particular, 
the energies p 2 + q 2 = -|n 2 + 2/ 2 are closer to the base 
energy -^n 2 than to the base energy of the next represen- 
tation, e n+ i » |n 2 + n. Of course, if the perturbation is 
symmetric under rotation, the next base energies coupled 
are e„ ±4 sa e„ ± 4n. 

There are, however, other unperturbed states with 
p 2 + q 2 w e„. For example, 7 2 + 49 2 = €70- However, the 
matrix elements of a smoothly perturbed Hamiltonian, 
7~ipq,p'q'-, in the unperturbed basis, are small if \p— p'\ 
or \q — q'\ is large. Because the perturbation due to a 
uniform field has a singular third derivative, these per- 
turbations drop off as a power law, and this relatively 
long range effect in momentum space is the source of the 
diffraction corrections. 

In Fig. 11 we show the magnitude of the amplitudes of 
the unperturbed states combining to make the state n = 
62, B = 25, m = 14. The area of a circle is proportional 
to the square of the amplitude. Above the main diagonal 
we show the theoretical result. All amplitudes lie on the 
line 31 — I, 31 4- I. Below, we show the result for the 
numerical wavefunction. The circle p 2 + q 2 = 2 • 31 2 is 
also shown. Our theory is starting to need corrections 
for m this large. 




FIG. 10. Streamlines for a sequence of states, 

m — 6,8, 10, 12, n = 62, B = 25 counterclockwise from the 
lower right. 

Thus, an interpretation of our method which yields the 
(1,1) resonance states of Eq. ([l2|) is that we effectively 
diagonalize the Hamiltonian in a basis restricted to the 
unperturbed states nearly 'degenerate' with e n and close 
to gn, hn. This is the case for the uniform field, and 



indeed, we achieve agreement between full numerical di- 
agonalization, diagonalization restricted to 'degenerate' 
states, and the procedure using the solution of the differ- 
ential equation Eq. (H). 

VII. AHARONOV-BOHM FLUX LINE 

The above approach can be generalized to deal with 
nonuniform flux configurations. To get the potential V 
associated with the (1, 1) resonance, all that is needed 
is to be able to calculate the flux contained within a 
(1,1) periodic orbit. We consider here the case of the 
Aharonov-Bohm flux line [ABFL]. Some further results 
are published elsewhere [ffq . 
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FIG. 11. The decomposition of the state of Fig. 6 in the 
5 = basis states. These are sin7rp(:E + |) sin irq(y + i) 
symmetrized and normalized. 

While much of the above discussion can be carried over 
to the ABFL, in the ideal case of a zero radius line there 
are strong diffraction effects which limit the applicability 
of our theory. We therefore begin by defining the alter- 
native problem of a finite size flux line or tube. This is 
more realistic if actual experiments are contemplated. 

Let p give the linear scale of the flux tube. We may 
think of the flux as uniform inside a tube of this radius, or 
as being distributed in some way, say as a gaussian, with 
p giving the scale of the distribution. [We actually used 
a square tube of side p. Even more accurately, we used 
four quarter strength tubes symmetrically located which 
allows a symmetry reduction in the numerics. In the ap- 
proximation of our theory, this gives essentially the same 
result as a single circular tube. The field inside a single 
flux tube is Bq — (fif-ip 2 and <j> is the total flux.] The 
typical angular deflection suffered by a classical particle 
traversing this field is 

66 ~ (4>/<fo)/kp. (23) 

To avoid diffraction we require 59 to be small. This can 
be achieved, of course, if kp is large, but that is not 
necessary. In the numerical work shown, we take <j)/<j)Q = 
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0.1, and p = 0.01, while k > 140. An alternative and 
equivalent condition is to insist that, on the appropriately 
defined average Q, the terms in the Hamiltonian satisfy 
({eA/c) 2 /2m) « (ep • A/mc). 

The results depend on where the flux line is located. 
We consider first the case that it is located at x = 0, 
y = — \ + a, where < a < \. Again, we consider states 
related to the (1,1) resonance. This leads as before to 
Eq. (|8|) but now we find a potential 



Vab(x) = —Bk/£; x G [—a, a], 

Vab(x) = +Bk/C; x G [-1 - a, -1 + a], 

Vab{x) — 0; x £ [-a, a] U [-1 - a, -1 + a] 



(24) 



for points in [—5,3], and Vab is extended periodically 
by V A b(x + 2) = Vab(x). For the flux line, B = 2t:<\>/^. 

Actually, Eq. ( |24| ) is for the ideal flux line. The ideal 
case has step function jumps in the potential which are 
smoothed out at finite p. In other words, the sharp jump 
at x = a, is replaced by a smooth rise beginning at x = 
a — p and ending at x — a + p. The exact shape of the 
rise depends on the distribution of flux in the line. As 
long as the transverse wavelengths in the solution of Eq. 
(^J) are long compared with p, the finiteness of p docs 
not play a significant role in our theory at this level of 
approximation. 

For the ideal, zero radius ABFL, there are significant 
deviations from this scenario. Indeed, most matrix el- 
ements of the ABFL perturbed Hamiltonian in the un- 
perturbed basis are infinite. However, it is a weak, log- 
arithmic infinity, and our theory seems to capture the 
main shape of the wavefunction, although at the rela- 
tively low energies for which numerical results are avail- 
able there are significant corrections. We consider these 
to be diffractive corrections, arising from a characteristic 
length shorter than the wavelength. 

It is clear that Bk/C is the important parameter in 
the UF case, while for the ABFL, both Bk/C and a are 
important. We begin with the case a = j, which has a 
'square well potential' of width \ near x = 0. For suf- 
ficiently large Bk/C, there will be 'tight binding' solu- 
tions approximately u m {x) = cos(m + l)irx/2a, \x\ < a, 
and zero elsewhere. This expression holds for sufficiently 
small even m, and for odd m the cosine is replaced by 
the sine. The energy E m » -Bk/C + n 2 (m + l) 2 /4a 2 . 

Fig. 12 is for the ABFL case, with n = 86, m = 0, 
r = and a = j. Fig. 12(a) shows u (x) and u (— x— 1), 
its extension into x < — i, reflected. For these pa- 
rameters, uq is not extremely localized, and extends 
significantly outside [— h,h]- The remai ning plots give 
\l/) m \ = \dV/dn\. Fig. 12(b) plots Eq. @ and 12(d) 
is from diagonalization in the limited basis of the 'de- 
generate' states. Clearly, these two approximations are 
nearly the same. Fig. 12(c) and 12(e) are obtained by 
numerical diagonalization in the complete unperturbed 
basis, for the finite size flux tube and the ideal flux line 
respectively. Clearly, for these parameters, there is sig- 
nificant diffraction from the flux tube. It is somewhat 
less than for the zero radius line but the two patterns 
have some resemblance. Diffraction evidently modifies 



the interference between different parts of the wavefunc- 
tion in an irregular way. The overall shape of \d^ /dn\ 
is well predicted by the theoretical uq. To give a sense 
of the irregularity of diffraction effects on the wave func- 
tion, we show in Fig. 12(f) the function for n = 70, a 
somewhat lower energy. In this case, according to theory, 
the uq is essentially the same, and the interference fringes 
have a somewhat longer wavelength, on the average. The 
diffractive effects are quite different in detail, however. 

It is clear from this result that the ABFL localizes the 
particle, and to do so, it must of course exert a force. 
That a flux line can exert a transverse force is now well 
established @. 

In Fig. 13 we show currents from a state with a = i, 
i.e. the much studied case with the flux line at the center 
of the square 0. Again, m = 0, r = 0, and n = 82. The 
upper part of the figure shows numerical streamlines for 
a corner of the square. Although the state is spatially 
not well localized, it does have a strongly 'paramagnetic' 
current structure. By paramagnetic we mean dE m /d<j) < 
0. 

The lower part shows the theoretical current den- 
sity j y (x) for y — 0, x G [— 4,0], that is, the 
current density along the upper edge of the up- 
per figure. A simple approximate formula for j y is 

j y (x) oc — (cosn7rx/2) 2 u m ,(— x — ^) 2 — u m (x — 

which gives double zeroes of the current at equally spaced 
x = (21 + l)/n. In this case j y is negative and the fac- 
tor depending on u m has no zeroes except at x = 0, for 
m = 0. 

The maxima are at j y = 0; j y < 0. Theory and nu- 
merics closely agree on the period and shape of the os- 
cillations, but the lower envelope of the two differ. We 
mark with horizontal bars the lower envelope of the nu- 
merical calculation, which is considerably more irregular 
than given by our theory. This structure of zeroes of j x is 
a consequence of a symmetry, namely y <-> —y, together 
with complex conjugation. This means that ^(x,0) is 
real and therefore gcnerically will have zeroes as a func- 
tion of x. 



VIII. SUMMARY 

A. Extensions of the results 

We have shown how to classify and find eigenstates of 
a charged particle in a square billiard subjected to a mag- 
netic flux which is classically weak. It is not necessarily 
weak quantally, however, and large remarkable changes 
in the wave functions are found which sometimes signifi- 
cantly localize the wave function. We used two basic flux 
configurations: a uniform field and an Aharonov-Bohm 
flux line. The placement of the ABFL is important. Of 
lesser importance is any finite radius to the line, at least 
in the range of parameters we use. We exhibited some 
of the wavefunctions in several forms and sequences, con- 
centrating on wavefunctions connected with the (1,1) res- 
onance, that is, connected with periodic orbits of the flux 
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free square whose velocities make angles of 45° with the 
coordinate axes. 

These results can be readily generalized to higher order 
resonances. For a given energy the effects diminish quite 
rapidly as the order increases, but in principle, at suffi- 
ciently high energy any given resonance can show strong 
magnetic effects. We give a few results in Appendix B. 
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FIG. 12. ABFL case, a =\,<j)= O.l^o, P = 0.01, n = 86, 
m = 0. a) Upper curve, uq(x), lower curve, u( — 1 — x). This 
is repeated in the other pictures. The remaining plots are of 
\d^/dn\ . b) is the theoretical result from Eq. Q13) , d) is from 
diagonalization in the 'degenerate' basis of Section V (E), and 
c) is from full numerical diagonalization. The interference 
oscillations scale as it(— 1 — x). e) and f) are numerical results 
for the ideal single ABFL (p = 0), with n — 82, and n — 70. 
The additional detailed structure is attributed to diffraction, 
but the overall scale is given well by our theory. 

One can extend these results to integrable systems 
other than the square billiard, and to other flux configu- 
rations. There are other billiard shapes, such as the rect- 
angle, certain triangles, the circle and the ellipse. One 
may also study 'soft billiards', e.g. a confining potential 
of the form U(x, y) — U\{x) + U2{y) where Ui is some sort 
of anharmonic potential. Soft billiards are more accurate 
representations of mesoscopic systems than are hard wall 
billiards, but the additional theoretical effort they require 
is not usually made. A practical difficulty, although not 
one of principle, is that it may be necessary to resort to 
action-angle variables, and it may be tedious to find the 
periodic orbits of the flux-free system, if some tiling trick 
cannot be used. 

In any case, the first steps of the program are to choose 
a convenient surface of section, choose the resonance of 
interest, find the AB phase &(x,x'), and the effective 
potential V(x). Much insight can be gained at this level. 

A much larger class of systems solvable by this tech- 
nique are nearly integrable systems subjected to a flux. 
For example, one could start from a nearly square trape- 
zoid. Then one would have a double perturbation of the 
square, one from the flux, the other from the change of 
shape. Each of these perturbations can have big quan- 
tum effects, and the combination of the two can be quite 
different from each one separately, especially since one 



breaks time reversal invariance and the other doesn't. 
Again, finding the effective potential is key to under- 
standing the qualitative results. 



B. Possible experiments 

There are possible experiments and even interesting 
devices which might be made, if time reversal can be 
broken. Having unusual wavefunctions suggests some of 
the possibilities. Probes of the system will depend on 
whether the wavefunction is large or small at the position 
of the probe. 




-0.5 -0,4-0.3 -0.2 -0,1 0.0 
x 

FIG. 13. Ideal zero radius flux line at the square center. 
Upper figure, numerical current streamlines for a quarter of 
the square. Lower figure, theoretical current density j y (x,0). 
Dashes indicate the numerical minima. 

Mesoscopic systems are of great current interest. In 
this case, electrons are the particles and there is a mag- 
netic flux. The system of Levy et al. || consisting of 
isolated metallic 'squares' is a case of this type, but one 
can imagine leads weakly connected to a metallic square 
whose coupling depends on the shape of the wavefunction 
at the Fermi level and the position of the lead. 

Another type of mesoscopic system is formed of surface 
electrons in 'corrals' on a metallic surface f2C|| . The wave 
functions of such electrons can be probed with atomic ac- 
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curacy with a scanning tunnelling microscope. Achieving 
an appropriate parameter range will be difficult, but per- 
haps not impossible. The corrals are leaky and do not 
really confine the electrons to their interior, but in many 
cases that idea seems to work, at least qualitatively. 

Another kind of system is the shallow square or rect- 
angular container of liquid, which is vibrated to set up 
standing waves. Instead of a flux, the system can be 
rotated with a uniform angular velocity. This differs 
from the uniform magnetic field in that the A 2 is ab- 
sent from the Hamiltonian, but since under our assump- 
tions that term can be neglected anyway, the same kind 
of results are expected. One can also introduce a non- 
ideal ABFL. This is done by making a small hole in the 
tank and allowing the water to flow out with a certain 
vorticity. Experiments on such a system have been per- 
formed although scattering rather than eigenstates 
was studied. Such an experiment would certainly have 
pedagogical value, and if the experiment is done, it may 
in fact be the first detailed observation of a nontrivial 
persistent current state. 

Still another system is the thin square microwave cav- 
ity. Here the 'quantum' waves are microwaves. There 
is no quantity directly equivalent to the magnetic field. 
However, one can replace part of the boundary by a fer- 
rite strip p2| , say between [—a, a] on one side of the cav- 
ity. The phase shift of the microwave upon reflection 
from the ferrite depends on the direction of magnetiza- 
tion of the ferrite and the direction of incidence of the 
microwave. This gives as an effective potential exactly 
that of Eq. (|2~4|). There is then a localized eigenstate cir- 
culating the cavity in one direction, but not the other. A 
similar situation from the point of view of diffraction and 
effective potential, but maintaining time reversal invari- 
ance, is the 'step' billiard p3| . This is a square billiard 
of side L say, with one side moved down by 5L = eL 
for \x\ < a. The perturbed (1,0) resonance states 'see' a 
square well effective potential to first approximation. If 
k5L is small, but kyLSL > 1, our theory works, gives 
nontrivial localization, and is similar to the flux line with 
small (f>/4>o. The range kSL — n, i.e. 6L a half wave- 
length, gives no phase shift for normal incidence and is 
similar to an ABFL with ip/^o — 1- Experiments on this 
system could be carried out in several contexts. 



C. Conclusions 

We have solved a characteristic example of quantum 
states in a weakly perturbed integrable system. The 
states and their energies can be classified and found to 
good approximation. The wavefunctions are quite non- 
trivial and interesting, as compared with the states of 
the integrable system. The states of a hard chaotic sys- 
tem are also not individually interesting, except for some 
weak 'scars', requiring resort to statistical studies and 
averages over many states in such systems. 

Our technique applies to a very large class of sys- 
tems, which includes experimental systems. The partic- 
ular case we emphasize breaks time reversal invariance, 



and the eigenstates have persistent currents. As far as 
we know, these are the first published nontrivial exam- 
ples such states. 
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APPENDIX A: TWO DIMENSIONAL 
WAVEFUNCTIONS 

In this Appendix we find the two dimensional wave- 
function. If a solution to = Tip is known, the wave- 
function in the square can be found by evaluating the 
integral over the surface of section 0| . 



*(r) 



G(r, x';E)ip(x')dx'. 



(Al) 



In this case 

G(r,x';E)= £ \ 



cl. t 

x 



d 2 S 



1/2 



27rfc dx'di± 

r 

exp |i<L?(r, x ; E) — i—v 



(A2) 



in our units. Here the sum is over the (unperturbed) 
classical trajectories going from a point x' on the sur- 
face of section to a point r inside the square, S(r,x';E) 
is the reduced action for this trajectory, £± is the direc- 
tion perpendicular to the trajectory at point r, and v 
is the Maslov index. In a single square scheme we take 
the surface of section to be the lower side of the square 
-1/2 <x'< 1/2, y' = -1/2. Then ip(x') is a linear com- 
bination of functions ipi and ipn described in the text 
such that i/>(±l/2) = (there should be no flow through 
the ends). It follows then that u m (x' + 2) = (— l) r u m (x') 
and 



i>(x') = e lKX 'u m (x') + i r e- iKx u m {x' - 1) 



(A3) 



with r defined by Eq. (13). Here we have also used a zero 
field quantization condition n — im/2 and the property 
Um(x') = (—l) m u m {—x'). Note that in the surface of 
section picture r does not have a direct interpretation as 
a representation label. 



We evaluate the integral (Al) in the stationary phase 
approximation. As we shall see, only the 45° orbits sur- 
vive. For a given point r = (x, y) there are four such 
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orbits: two starting at x[ and two at x' 2 (Fig. 14). Each 
orbit gives one term in Eq. (|l2|). Note that the orbits 
with the positive (negative) ^-projection of the momen- 
tum at x' are generated by the first (second) term in Eq. 



(A3) 
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FIG. 14. A diagram showing how the full wavefunction 
ty(x,y) is constructed from the surface of section wavefunc- 
tion. Two points on the surface of section, x' x , x' 2 dominate in 
the sense of stationary phase, and there are two paths from 
each of these points to the point (x,y). 

Although the following calculation depends in which 
sector inside the square the point lies, the final result 
is independent of sector. So for definitiveness we as- 
sume y < x < —y, that is the point lies below the di- 
agonals y — x and y — —x. We have numbered the 
orbits by s = 0, . . . , 3 according to the terms in Eq. 
( |l2| ) they represent. For s = orbit, expressing the ac- 
tion in terms of the distance, S(r,x';E) = kL(r,x') = 
k\/ (x — x') 2 + (y + 1/2) 2 , we find the prefactor in Eq. 



done by the rule: s = trajectory arrives to r from the 
South- West, s = 1 from SE, s — 2 from NE, and s = 3 
from NW. 

It is worth noticing that Eqs. (A5) or ( |l2] ) are valid 
only up to the square root order of the small parameter 
e = B/k of the perturbation theory behind our work. 
If we wish to obtain the results valid to the first or- 
der of e we should (a) use a better approximation for 
Um(x') [fy], (b) add the vector potential term to the ac- 
tion S(r^x; E), and (c) find a correction to the stationary 
point due to u m {x'). 

Alternatively, the above calculation can be carried out 
in the repeated square scheme defined in Sec. III. The 
surface of section is then the line y = —1/2 (identified 
with y = 3/2). We can restrict ourselves to the tra- 
jectories with the positive x- and y-projections of the 
momentum in this extended manifold. Then only ipi(x') 
is needed to generate a complete two-dimensional wave- 
function. In order for i/jj to have a period 2, the condition 
u m (x! + 2) = (— l) r u m (x') must be satisfied. A point r 
in the original square will be represented by its images 
in four (up to a translation by Ax = 2) domains in the 
manifold. Each of these images is connected by one 45° 
orbit to the surface of section. These four orbits generate 
the four terms s = 0, . . . , 3 in Eq. ( |l2] ) as follows: s = is 
produced by the point r in the original domain, s = 1, 3, 2 
by its reflection about the right side of the square, upper 
side of the square, and the composition of both respec- 
tively. There are no Maslov phases in Eq. (A.2) in this 
case. However the terms s = 1, 3, which are obtained by 
the odd number of reflections, should be summed with 
an additional minus sign. Indeed, the Jacobian of the co- 
ordinate transformation from the original square to the 
extended manifold is singular on the boundary. So, when 
the wavefunction on the manifold is folded back to the 
physical domain a phase difference n will be accumulated 
between each pair of domains related by one reflection. 
Of course, this phase is analogous to the Maslov phase 
for an impenetrable wall. 
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(A4) 



The stationary point in the integral (Al) is determined 
by t he exponents in Eq. ( A2) and in the first term of Eq. 
(|A|). It is x[ = x- y- 1/2, as in a 45° orbit. The s = 
wavefunction is then 



_ J.K(x+y) 



1 

x y 

y 2 



(A5) 



dropping a constant factor K~ 1 / 2 e 4 ( K / 2 ~ 7r//4 V2. The 
terms with s = 1,2,3 can be obtained in the same man- 
ner with an appropriate addition of the Maslov phase for 
each bounce. 

If r is located in a different sector of the square, the 
results remain the same if proper account is made of the 
Maslov phases. The labeling of the orbits is invariant if 



APPENDIX B: HIGHER RESONANCES 



For completeness, we give some results for higher res- 
onances, which can be labelled with relatively prime in- 
tegers, p, q. These resonances are closely related to states 
of a rectangle of side 1/p, 1/q. Obviously, for the square, 
the energies for p, q are the same as those for q, p. The 
corresponding states are also the same after a 90° rota- 
tion. However, for \p — q\ / \p + q\ not too small, the p, q 
states and the q,p states are different and are nearly un- 
coupled and there is a nearly doubly degenerate set of 
energy levels. The splitting of the exact levels, which are 
combinations of pq and qp that are eigenstates of the ro- 
tation operator 1Z, can be estimated by using an analog 
of 'chaos assisted tunnelling'. 

The classical periodic orbits make q bounces from the 
x sides and p bounces from the y sides. These resonant 
orbits all have the same length C pq = 2^/p 2 + q 2 . The 
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maximal directed area enclosed by an orbit is ±l/2pq. 
The energy is k\ m where 



2irn 

7* 



1+^ 

9 . 



(Bl) 



and n is an integer. The first term is much larger than the 
second. Again, we have the approximate quantization of 
n wavelengths in an orbit length. 

The potential V pq (x) depends on q being associated 

with the x-direction, thus the notation Em . This 'trans- 
verse' contribution to the energy is of course symmetric 
in p and q, that is, q 2 E$ — p 2 E$ . The potential is 



V Pq {x) 





3 






1" 


V 


q ^x-\ 
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(B2) 



where V is given by Eq. (g). V pq has period 2/q, rather 
than 2. The boundary condition on the eigenstates is 
u(x — 2) = e %l3 u(x) where (3 = 27rfrac \pn/ (p 2 + q 2 )] 
and frac indicates the fractional part. The potential can 
be regarded as weak if kB << pqC pq . If this condition is 
satisfied the pq resonance can be ignored. 

Similar but somewhat more complex results can be 
obtained for rectangles and certain triangles. 
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